#!/usr/bin/perl -w
use strict;
use FindBin;
use lib ("$FindBin::Bin/..", "/net/cpp-group/Leo/bin");
use Getopt::Long;
use parse_bl2seq;
use db_parameters;
use File::Temp qw/ tempfile/;
#88888888888888888888888888888888888888888888888888888888888888888888888888888888

#   Usage

#88888888888888888888888888888888888888888888888888888888888888888888888888888888
my $usage = <<'USAGE';

USAGE:

    get_aligned_sequences.pl
                            [--a_species "human"]
                            [--b_species "mouse"]
                            [--cdna_instead_of_protein]
                            [--help]
                            IDs.data


        This script takes a file of tab delimited gene_ids via
        STDIN or on the command line and prints the aligned clustal sequence
        pairs to STDOUT
USAGE

#88888888888888888888888888888888888888888888888888888888888888888888888888888888888888888

#   Get options

#88888888888888888888888888888888888888888888888888888888888888888888888888888888888888888
# mandatory


# optional parameters
my $help = undef;
my $a_species = 'human';
my $b_species = 'mouse';
my $cdna_mode = undef;
GetOptions('help' => \$help,
		   'a_species:s' => \$a_species,
		   'b_species:s' => \$b_species,
		   'cdna_instead_of_protein' => \$cdna_mode,
		   );

die $usage if ($help);
die $usage unless $a_species;
die $usage unless $b_species;
die $usage unless @ARGV;





my ($fh, $filename);

# delete files automatically when we exit or die
END
{
	unlink($filename) if defined $filename;
}




eval {
	($fh, $filename) = tempfile( 'testXXXXXX', SUFFIX => '.dat'); #, UNLINK => 1
	print $fh <<'PL/SQLCMD';
		CREATE TEMP TABLE t_1 (gene_id1 TEXT, gene_id2 TEXT);
		COPY t_1 from STDIN;
PL/SQLCMD
	my $cnt_sequences;
	while (<>)
	{
		chomp;
		next unless /\t/;
		$cnt_sequences++;
		print $fh $_, "\n";
	}
	my $table_name = 'orthologs.'.$a_species . '_' . $b_species . '_sequences';
	print $fh <<"PL/SQLCMD";
\\.
		SELECT gene_id1, gene_id2, sequence1, sequence2, cdna1, cdna2, gaps1, gaps2, align_start, align_len
		FROM $table_name
		NATURAL JOIN t_1;
PL/SQLCMD
	close $fh;
	my $psql_cmd = "psql -t -A -q -F \"\t\" -h $db_host -d $db_name -U $db_user -f $filename |";
	open (FH_PSQL, $psql_cmd)
		or die "Error:\n\tCould not open pipe to psql\n";


	use constant ENSG1 => 0;
	use constant ENSG2 => 1;
	use constant PROT1 => 2;
	use constant PROT2 => 3;
	use constant CDNA1 => 4;
	use constant CDNA2 => 5;
	use constant GAPS1 => 6;
	use constant GAPS2 => 7;
	use constant START => 8;
	use constant LENTH => 9;
	print "CLUSTAL W\n";
	while (<FH_PSQL>)
	{
		chomp;
		my @data = split /\t/;
		my ($seq1, $seq2);
		if ($cdna_mode)
		{
			$seq1 = substr(expand_cdna_with_alignment_str($data[CDNA1], "-", $data[GAPS1]),
							$data[START] * 3, $data[LENTH] * 3);
			$seq2 = substr(expand_cdna_with_alignment_str($data[CDNA2], "-", $data[GAPS2]),
							$data[START] * 3, $data[LENTH] * 3);
		}
		else
		{
			$seq1 = substr(expand_prot_with_alignment_str($data[PROT1], "-", $data[GAPS1]),
							$data[START], $data[LENTH]);
			$seq2 = substr(expand_prot_with_alignment_str($data[PROT2], "-", $data[GAPS2]),
							$data[START], $data[LENTH]);
		}
		print $data[ENSG1], "\t", $seq1, "\n",
			  $data[ENSG2], "\t", $seq2, "\n";
	}
};
if ($@) {
	print STDERR "unlink($filename)\n";
	unlink($filename) if defined $filename;
	$filename = undef;
}

